import numpy as np
import csv
import random
from datetime import datetime
from time import strftime
from scipy.stats import scoreatpercentile
from dbfpy import dbf

'''
revision of test24
to calculate the false nagetive, false positive

'''

def getCurTime():
    """
    get current time
    Return value of the date string format(%Y-%m-%d %H:%M:%S)
    """
    format='%Y-%m-%d %H:%M:%S'
    sdate = None
    cdate = datetime.now()
    try:
        sdate = cdate.strftime(format)
    except:
        raise ValueError
    return sdate

def build_data_list(inputCSV):
    sKey = []
    fn = inputCSV
    ra = csv.DictReader(file(fn), dialect="excel")
    
    for record in ra:
        #print record[ra.fieldnames[0]], type(record[ra.fieldnames[-1]])
        for item in ra.fieldnames:
            temp = float(record[item])
            sKey.append(temp)
    sKey = np.array(sKey)
    sKey.shape=(-1,len(ra.fieldnames))
    return sKey

def build_region_list(inputCSV):
    temp =  build_data_list(inputCSV)
    temp_list = temp[np.argsort(temp[:,0]),:]
    #temp_idset = set(temp_list[:,1])
    temp_idDict = {}
    i = 0
    for item in temp_list:
        if str(int(item[1])) not in temp_idDict.keys():
            temp_idDict[str(int(item[1]))] = i
            i += 1
    #print i
    temp_idlist = []
    for item in temp_list[:,1]:
        temp_idlist.append(temp_idDict[str(int(item))])
    return temp_idlist

def cal_region_attri(region_id):
    dis_reg = np.unique(region_id)
    iLen = dis_reg.shape[0]
    temp_reg_attri = np.zeros((iLen, 3)) #[caner, pop, rate]

    # unit_attri = [cancer, pop, area, rate]

    i = 0
    for item in unit_attri:
        temp_reg_attri[int(region_id[i]),0] += item[0]
        temp_reg_attri[int(region_id[i]),1] += item[1]
        #temp_reg_attri[int(region_id[i]),2] + = item[2]
        i = i + 1

    for item in temp_reg_attri:
        item[2] = item[0]/item[1]

    return temp_reg_attri[:,2]

def cal_unit_rate(inputCSV):           
    region_id = build_region_list(inputCSV)
    reg_rate = cal_region_attri(region_id) # [rate]
    rate = []
    for item in region_id:
        rate.append(reg_rate[int(item)])
    return rate, region_id

def build_unit_pvalue(inputCSV):
    unit_pvalue = []
    fn = inputCSV
    ra = csv.DictReader(file(fn), dialect="excel")
    
    for record in ra:
        #print float(record[ra.fieldnames[-1]]), type(record[ra.fieldnames[-1]])
        unit_pvalue.append(float(record[ra.fieldnames[-1]]))
    return unit_pvalue

def build_satscanresult_dbf(inputDBF):
    sKey = np.array([])
    fn = inputDBF
    db = dbf.Dbf(fn)
    #recordLenth = 0
    for record in db:
        temp = float(record[db.fieldNames[2]])
        temp_id = int(float(record[db.fieldNames[0]]))
        sKey = np.append(sKey, [temp_id, temp])

    sKey.shape = (-1,2)
    return sKey

def cal_quantile(inputdata, q):
    temp = []
    #print len(inputdata[0])
    for i in range(0, len(inputdata[0])):
        #print scoreatpercentile(inputdata[:,0], 25)
        temp.append(scoreatpercentile(inputdata[:,i], int(q), limit = ()))
    return temp

def cal_riskare_attri():
    i = 0
    temp = [0, 0, 0, 0, 0, 0]
    for item in unit_attri:
        if i in high_risk_area_id:
            temp[0] += item[0]
            temp[1] += item[1]
        elif i in low_risk_area_id:
            temp[2] += item[0]
            temp[3] += item[1]
        i += 1
    temp[4] = temp[0] + temp[2]
    temp[5] = temp[1] + temp[3]
    temp = np.array(temp)
    temp.shape = (3,-1)
    # [[high_case, high_pop]
    # [low_case, low_pop]
    # [total_case, total_pop]]
    return temp

def cal_TP():
    pvalueSet = set(unit_attri[:,2])
    pvalueList = []
    for item in pvalueSet:
        pvalueList.append(item)
    pvalueList.sort()
    #print pvalueSet

    maxHighCaseMeasure = np.zeros(4)
    maxLowCaseMeasure = np.zeros(4)
    maxHighPopMeasure = np.zeros(4)
    maxLowPopMeasure = np.zeros(4)
    maxCaseMeasure = np.zeros(4)
    maxPopMeasure = np.zeros(4)

    perfect_p = np.zeros(4) # highCase, lowCase, highPop, lowPop
    #print unit_attri[:,2]
    #print pvalueList
    
    for p in pvalueList:
        tempHighCaseMeasure = np.zeros(4)
        tempLowCaseMeasure = np.zeros(4)
        tempHighPopMeasure = np.zeros(4)
        tempLowPopMeasure = np.zeros(4)
        tempCaseMeasure = np.zeros(4)
        tempPopMeasure = np.zeros(4)

        i = 0
        for item in unit_attri:
            if item[2] < p:
                if i in high_risk_area_id:
                    tempHighCaseMeasure[0] += item[0]
                    tempHighPopMeasure[0] += item[1]
                else:
                    tempHighCaseMeasure[1] += item[0]
                    tempHighPopMeasure[1] += item[1]
                if i in low_risk_area_id:
                    tempLowCaseMeasure[0] += item[0]
                    tempLowPopMeasure[0] += item[1]
                else:
                    tempLowCaseMeasure[1] += item[0]
                    tempLowPopMeasure[1] += item[1]
                if i in risk_area_id:
                    tempCaseMeasure[0] += item[0]
                    tempPopMeasure[0] += item[1]
                else:
                    tempCaseMeasure[1] += item[0]
                    tempPopMeasure[1] += item[1]
            i += 1

        #print riskarea_attri
        tempHighCaseMeasure[2] = riskarea_attri[0,0] - tempHighCaseMeasure[0]
        tempHighPopMeasure[2] = riskarea_attri[0,1] - tempHighPopMeasure[0]
        tempLowCaseMeasure[2] = riskarea_attri[1,0] - tempLowCaseMeasure[0]
        tempLowPopMeasure[2] = riskarea_attri[1,1] - tempLowPopMeasure[0]
        tempCaseMeasure[2] = riskarea_attri[2,0] - tempCaseMeasure[0]
        tempPopMeasure[2] = riskarea_attri[2,1] - tempPopMeasure[0]

        tempHighPopMeasure[3] = (tempHighPopMeasure[0] + 0.0)/(tempHighPopMeasure[1] + tempHighPopMeasure[2])
        tempHighCaseMeasure[3] = (tempHighCaseMeasure[0] + 0.0)/(tempHighCaseMeasure[1] + tempHighCaseMeasure[2])
        
        tempLowPopMeasure[3] = (tempLowPopMeasure[0] + 0.0)/(tempLowPopMeasure[1] + tempLowPopMeasure[2])
        tempLowCaseMeasure[3] = (tempLowCaseMeasure[0] + 0.0)/(tempLowCaseMeasure[1] + tempLowCaseMeasure[2])
    
        tempPopMeasure[3] = (tempPopMeasure[0] + 0.0)/(tempPopMeasure[1] + tempPopMeasure[2])
        tempCaseMeasure[3] = (tempCaseMeasure[0] + 0.0)/(tempCaseMeasure[1] + tempCaseMeasure[2])

        #print p, tempHighPopMeasure[3]
        
        #perfect_p: highCase, lowCase, highPop, lowPop

        if tempHighCaseMeasure[3] > maxHighCaseMeasure[3]:
            maxHighCaseMeasure = tempHighCaseMeasure
            perfect_p[0] = p
        if tempLowCaseMeasure[3] > maxLowCaseMeasure[3]:
            maxLowCaseMeasure = tempLowCaseMeasure
            perfect_p[1] = p
        if tempHighPopMeasure[3] > maxHighPopMeasure[3]:
            maxHighPopMeasure = tempHighPopMeasure
            perfect_p[2] = p
        if tempLowPopMeasure[3] > maxLowPopMeasure[3]:
            maxLowPopMeasure = tempLowPopMeasure
            perfect_p[3] = p
            
        if tempCaseMeasure[3] > maxCaseMeasure[3]:
            maxCaseMeasure = tempPopMeasure
            #perfect_p[1] = p
        if tempPopMeasure[3] > maxPopMeasure[3]:
            maxPopMeasure = tempPopMeasure
            #perfect_p[2] = p

    maxCaseMeasure = np.append(maxCaseMeasure, maxHighCaseMeasure)
    maxCaseMeasure = np.append(maxCaseMeasure, maxLowCaseMeasure)

    maxPopMeasure = np.append(maxPopMeasure, maxHighPopMeasure)
    maxPopMeasure = np.append(maxPopMeasure, maxLowPopMeasure)

    return maxCaseMeasure, maxPopMeasure, perfect_p

#--------------------------------------------------------------------------
#MAIN

if __name__ == "__main__":
    print "begin at " + getCurTime()

    H1 = [8,16,844,915,919,921,923,924]
    L2 = [5,103,106,513,517,518,520,531,534,535,536,541]
    H3 = [63,265,267,268,333,336,337,339,340,342,343,348]
    H4 = [13,174,178,198,886,887,888,889,890]
    L5 = [146,171,182,810,811,814,815,864,867]
    L6 = [20,133,692,694,695,696,698,702,705]
    H7 = [69,70,87,88,369,370,372,442,443]
    high_risk_area_id = H1 + H3 + H4 + H7
    low_risk_area_id = L2 + L5 + L6
    risk_area_id = high_risk_area_id + low_risk_area_id
    
    #unitCSV = "C:/TP1000_1m.csv"
    unitCSV = 'C:/_DATA/CancerData/test/Jun08/TP1000_1m_16_04.csv'
    dataMatrix = build_data_list(unitCSV)  # [id, pop, cancer1, cancer2, cancer3]
    unit_attri = np.zeros((len(dataMatrix),4))
    
    # unit_attri [cancer, pop, rate, id]
    unit_attri[:,1] = dataMatrix[:,1]
    #unit_attri[:,3] = np.arange(0,1000)
    perfect_pvalue = []

    j = 0
    while j < 1:
        #filePath = 'C:/_DATA/CancerData/test/Jan15/Thousand/10000/kernel_LLR/'
        #filePath = 'C:/Documents and Settings/wang322/Desktop/'
        filePath = 'C:/_DATA/CancerData/test/Jun08/satscan/highlow/'

        # returnID used in function find_rep_region
        returnID = 1  # 0: cancer, 1: pop, 2: count
        repeatTime = 1000
        while returnID < 2:
            case_measure = np.array([])
            pop_measure = np.array([])

            for repeat in range(0, repeatTime):
                #repeat += 999
                print repeat
                #regionCSV = filePath + 'LLR_EBS_high_' + str(repeat) +"_pvalue.csv"
                regionDBF = filePath + 'hl' + str(repeat) +".gis.dbf"
                #regionCSV = filePath + 'LLR_SSD_high_pvalue.csv'
                #regionCSV = filePath + "temp.csv"
                unit_attri[:,0] = dataMatrix[:,repeat + 2] # [cancer, pop, pvalue]
                unit_attri[:,2] = np.ones(len(unit_attri))
                temp_pvalue = build_satscanresult_dbf(regionDBF)
                #print temp_pvalue
                for item in temp_pvalue:
                    unit_attri[int(item[0]),2] = item[1]
                #print unit_attri[:,2]
                riskarea_attri = cal_riskare_attri()  # [total_cancer, total_pop, total_area]
                #print unit_attri[:,2]
                
                temp_case_measure, temp_pop_measure, temp_p = cal_TP()
                '''
                
                temp_case_measure = np.append(temp_case_measure, temp_case_measure[0] + temp_case_measure[4])
                temp_case_measure = np.append(temp_case_measure, temp_case_measure[1] + temp_case_measure[5])
                temp_case_measure = np.append(temp_case_measure, temp_case_measure[2] + temp_case_measure[6])
                temp_case_measure = np.append(temp_case_measure, temp_case_measure[8]/(temp_case_measure[9] + temp_case_measure[10]))
                
                temp_pop_measure = np.append(temp_pop_measure, temp_pop_measure[0] + temp_pop_measure[4])
                temp_pop_measure = np.append(temp_pop_measure, temp_pop_measure[1] + temp_pop_measure[5])
                temp_pop_measure = np.append(temp_pop_measure, temp_pop_measure[2] + temp_pop_measure[6])
                temp_pop_measure = np.append(temp_pop_measure, temp_pop_measure[8]/(temp_pop_measure[9] + temp_pop_measure[10]))
                '''
                
                case_measure = np.append(case_measure, temp_case_measure)
                pop_measure = np.append(pop_measure, temp_pop_measure)
                
                perfect_pvalue.append(temp_p)
            #print 'here'
            case_measure.shape = (repeatTime, -1)
            
            
            temp_mean = case_measure.mean(axis=0)
            temp_1Q = cal_quantile(case_measure, 25)
            temp_median = np.median(case_measure, axis=0)
            temp_3Q = cal_quantile(case_measure, 75)
            #temp_std = case_measure.std(axis=0)
            
            case_measure = np.append(case_measure, temp_mean)
            case_measure = np.append(case_measure, temp_1Q)
            case_measure = np.append(case_measure, temp_median)
            case_measure = np.append(case_measure, temp_3Q)
            #case_measure = np.append(case_measure, temp_std)
            case_measure.shape = (repeatTime + 4, -1)
            
            pop_measure.shape = (repeatTime, -1)
            #scoreatpercentile(pop_measure[:0], 25, limit = ())
    
            temp_mean = pop_measure.mean(axis=0)
            temp_1Q = cal_quantile(pop_measure, 25)
            temp_median = np.median(pop_measure, axis=0)
            temp_3Q = cal_quantile(pop_measure, 75)
            #temp_std = pop_measure.std(axis=0)
            pop_measure = np.append(pop_measure, temp_mean)
            pop_measure = np.append(pop_measure, temp_1Q)
            pop_measure = np.append(pop_measure, temp_median)
            pop_measure = np.append(pop_measure, temp_3Q)
            pop_measure.shape = (repeatTime + 4, -1)
            #print pop_measure[:,3]
            
            #print pop_measure[:,1].max(), count_measure[:,1].max(), pop_measure[:,6].max(), count_measure[:,6].max(), pop_measure[:,1].min(),count_measure[:,1].min(), pop_measure[:,6].min(), count_measure[:,6].min()
            #fileLoc = filePath + 'high_case_measure_1.csv'
            #np.savetxt(fileLoc, case_measure, delimiter=',', fmt = '%10.5f')
            fileLoc = filePath + 'pop_measure_1.csv'
            np.savetxt(fileLoc, pop_measure, delimiter=',', fmt = '%10.5f')
            fileLoc = filePath + 'pvalue_1.csv'
            np.savetxt(fileLoc, perfect_pvalue, delimiter=',', fmt = '%10.5f')
            #print pop_measure
            #print pop_measure[:,-1]
            returnID = returnID + 1
        j = j + 1

    print "end at " + getCurTime()
    print "========================================================================"  

           
